Ocean acidification and boring sponge on Crassostrea virginica GE
Last run on 29 June 2022
Expoloring data
Filtering counts
After viewing the size factors, it appears that we have 5 samples with really low values that should be removed:
- 19355
- 19372
- 19382
- 19496
- 19505
Without filtering low reads, we had a total of 38828 counts. After filtering out the low counts (those with a base mean less than 3), we now have 24384 counts remaining for all C. virginica samples.
PCA plots of Global Gene Expression
# going to use the filter_counts_out object created above because it is filtered for low reads and has outliers removed
## create the DESeq object
dds <- DESeqDataSetFromMatrix(countData = filter_counts_out,
colData = expDesign,
design = ~ treat)
dds <- DESeq(dds) # differential expression analysis on gamma-poisson distribution
vsd <-varianceStabilizingTransformation(dds, blind = TRUE) # quickly estimate dispersion trend and apply a variance stabilizing transformation
## Save dds and vsd data
save(dds, vsd, file = "Data/transformed_counts.RData")## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 1500
##
## adonis2(formula = pcaData[6:31] ~ pcaData$infect * pcaData$pCO2, permutations = 1500, method = "eu")
## Df SumOfSqs R2 F Pr(>F)
## pcaData$infect 1 1118.3 0.04736 1.1893 0.1799
## pcaData$pCO2 1 945.3 0.04003 1.0053 0.3817
## pcaData$infect:pcaData$pCO2 1 864.0 0.03659 0.9188 0.5623
## Residual 22 20686.8 0.87603
## Total 25 23614.3 1.00000
Differential gene expression
DEGs
## log2 fold change (MLE): treat S 400 vs N 2800
## Wald test p-value: treat S 400 vs N 2800
## DataFrame with 24384 rows and 6 columns
## baseMean log2FoldChange lfcSE
## <numeric> <numeric> <numeric>
## LOC111126949 3.45225877157332 -0.942224682074193 0.817783710717942
## LOC111110729 0.485466216780156 -0.512448068281346 1.46364594751556
## LOC111120752 1.69362126099025 -0.282706315518707 0.84681189235072
## LOC111113860 5.93338041130341 -0.0448270626067275 0.601266321422372
## LOC111109550 0.331392203994242 -0.986053856261322 2.48365300058619
## ... ... ... ...
## LOC111117112 0.646302530637881 1.67887461151877 1.80337428132596
## LOC111117113 0.14363842589226 0.768137359775011 3.83879218552636
## LOC111117117 0.0778207961311771 0.287240683720725 4.24705363394974
## LOC111116908 0.109745657738521 -0.0817811312489327 4.24579857131686
## LOC111117715 0.621026843813011 -0.445504180182541 1.6799677252913
## stat pvalue padj
## <numeric> <numeric> <numeric>
## LOC111126949 -1.15216856208471 0.249251813595227 0.999651667882862
## LOC111110729 -0.35011750563802 0.726250513749595 0.999651667882862
## LOC111120752 -0.333847833352841 0.738494386321806 0.999651667882862
## LOC111113860 -0.074554421243291 0.940569239720562 0.999651667882862
## LOC111109550 -0.397017560838247 0.691354510914985 0.999651667882862
## ... ... ... ...
## LOC111117112 0.930962933709109 0.351872738132339 0.999651667882862
## LOC111117113 0.20009870882596 0.841403383153489 0.999651667882862
## LOC111117117 0.0676329306097302 0.946077840571377 0.999651667882862
## LOC111116908 -0.019261660645283 0.984632368623089 0.999651667882862
## LOC111117715 -0.265186154159772 0.790866060568095 0.999651667882862##
## out of 24380 with nonzero total read count
## adjusted p-value < 0.1
## LFC > 0 (up) : 2, 0.0082%
## LFC < 0 (down) : 0, 0%
## outliers [1] : 36, 0.15%
## low counts [2] : 3, 0.012%
## (mean count < 0)
## [1] see 'cooksCutoff' argument of ?results
## [2] see 'independentFiltering' argument of ?resultssummary(res_con_sponge)##
## out of 24380 with nonzero total read count
## adjusted p-value < 0.05
## LFC > 0 (up) : 0, 0%
## LFC < 0 (down) : 0, 0%
## outliers [1] : 36, 0.15%
## low counts [2] : 3, 0.012%
## (mean count < 0)
## [1] see 'cooksCutoff' argument of ?results
## [2] see 'independentFiltering' argument of ?resultssummary(res_con_pco2)##
## out of 24380 with nonzero total read count
## adjusted p-value < 0.05
## LFC > 0 (up) : 0, 0%
## LFC < 0 (down) : 0, 0%
## outliers [1] : 36, 0.15%
## low counts [2] : 3, 0.012%
## (mean count < 0)
## [1] see 'cooksCutoff' argument of ?results
## [2] see 'independentFiltering' argument of ?resultssummary(res_trt)##
## out of 24380 with nonzero total read count
## adjusted p-value < 0.05
## LFC > 0 (up) : 6, 0.025%
## LFC < 0 (down) : 5, 0.021%
## outliers [1] : 36, 0.15%
## low counts [2] : 11339, 47%
## (mean count < 2)
## [1] see 'cooksCutoff' argument of ?results
## [2] see 'independentFiltering' argument of ?results#summary(res_trt_pco2)
Venn of DEGs
It is just a demo, not very exciting :D
Session Information
Session information from the last full knit of Rmarkdown on 29 June 2022.
## R version 3.6.3 (2020-02-29)
## Platform: x86_64-apple-darwin15.6.0 (64-bit)
## Running under: macOS Catalina 10.15.7
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] grid parallel stats4 stats graphics grDevices utils
## [8] datasets methods base
##
## other attached packages:
## [1] ggvenn_0.1.9 adegenet_2.1.4
## [3] ade4_1.7-17 ggrepel_0.9.1
## [5] pdftools_3.0.1 ggpubr_0.4.0
## [7] data.table_1.14.0 vegan_2.5-7
## [9] lattice_0.20-44 permute_0.9-5
## [11] DESeq2_1.26.0 SummarizedExperiment_1.16.1
## [13] DelayedArray_0.12.3 BiocParallel_1.20.1
## [15] matrixStats_0.57.0 Biobase_2.46.0
## [17] GenomicRanges_1.38.0 GenomeInfoDb_1.22.1
## [19] IRanges_2.20.2 S4Vectors_0.24.4
## [21] BiocGenerics_0.32.0 plotly_4.9.4.1
## [23] arrayQualityMetrics_3.42.0 forcats_0.5.1
## [25] stringr_1.4.0 dplyr_1.0.7
## [27] purrr_0.3.4 readr_2.0.1
## [29] tidyr_1.1.3 tibble_3.1.3
## [31] ggplot2_3.3.5 tidyverse_1.3.1
## [33] knitr_1.33
##
## loaded via a namespace (and not attached):
## [1] utf8_1.2.2 tidyselect_1.1.1 RSQLite_2.2.7
## [4] AnnotationDbi_1.48.0 htmlwidgets_1.5.3 beadarray_2.36.1
## [7] munsell_0.5.0 units_0.7-2 codetools_0.2-18
## [10] preprocessCore_1.48.0 withr_2.4.2 colorspace_2.0-2
## [13] highr_0.9 rstudioapi_0.13 setRNG_2013.9-1
## [16] ggsignif_0.6.2 labeling_0.4.2 GenomeInfoDbData_1.2.2
## [19] hwriter_1.3.2 farver_2.1.0 bit64_4.0.5
## [22] LearnBayes_2.15.1 coda_0.19-4 vctrs_0.3.8
## [25] generics_0.1.0 xfun_0.30 qpdf_1.1
## [28] R6_2.5.1 illuminaio_0.28.0 gridSVG_1.7-2
## [31] locfit_1.5-9.4 bitops_1.0-7 cachem_1.0.5
## [34] assertthat_0.2.1 promises_1.2.0.1 scales_1.1.1
## [37] nnet_7.3-16 gtable_0.3.0 affy_1.64.0
## [40] rlang_0.4.11 genefilter_1.68.0 systemfonts_1.0.2
## [43] splines_3.6.3 rstatix_0.7.0 lazyeval_0.2.2
## [46] hexbin_1.28.2 broom_0.7.9 checkmate_2.0.0
## [49] BiocManager_1.30.16 yaml_2.3.5 reshape2_1.4.4
## [52] abind_1.4-5 modelr_0.1.8 crosstalk_1.1.1
## [55] backports_1.2.1 httpuv_1.6.2 Hmisc_4.5-0
## [58] tools_3.6.3 spData_0.3.10 affyio_1.56.0
## [61] ellipsis_0.3.2 raster_3.4-13 jquerylib_0.1.4
## [64] RColorBrewer_1.1-2 proxy_0.4-26 Rcpp_1.0.7
## [67] plyr_1.8.6 base64enc_0.1-3 zlibbioc_1.32.0
## [70] classInt_0.4-3 RCurl_1.98-1.4 deldir_0.2-10
## [73] rpart_4.1-15 openssl_1.4.4 haven_2.4.3
## [76] cluster_2.1.2 fs_1.5.0 magrittr_2.0.1
## [79] openxlsx_4.2.4 gmodels_2.18.1 reprex_2.0.1
## [82] hms_1.1.0 mime_0.11 evaluate_0.14
## [85] xtable_1.8-4 XML_3.99-0.3 rio_0.5.27
## [88] jpeg_0.1-9 gcrma_2.58.0 readxl_1.3.1
## [91] gridExtra_2.3 compiler_3.6.3 KernSmooth_2.23-20
## [94] crayon_1.4.1 htmltools_0.5.2 spdep_1.1-8
## [97] mgcv_1.8-36 later_1.3.0 tzdb_0.1.2
## [100] Formula_1.2-4 geneplotter_1.64.0 expm_0.999-6
## [103] lubridate_1.7.10 DBI_1.1.1 dbplyr_2.1.1
## [106] MASS_7.3-54 boot_1.3-28 sf_1.0-5
## [109] Matrix_1.3-4 car_3.0-11 cli_3.0.1
## [112] vsn_3.54.0 gdata_2.18.0 igraph_1.2.6
## [115] pkgconfig_2.0.3 sp_1.4-5 foreign_0.8-75
## [118] xml2_1.3.2 svglite_2.0.0 annotate_1.64.0
## [121] bslib_0.3.1 BeadDataPackR_1.38.0 affyPLM_1.62.0
## [124] XVector_0.26.0 rvest_1.0.1 digest_0.6.27
## [127] Biostrings_2.54.0 rmarkdown_2.13 base64_2.0
## [130] cellranger_1.1.0 htmlTable_2.2.1 curl_4.3.2
## [133] gtools_3.9.2 shiny_1.7.1 lifecycle_1.0.0
## [136] nlme_3.1-151 jsonlite_1.7.2 carData_3.0-4
## [139] seqinr_4.2-8 viridisLite_0.4.0 askpass_1.1
## [142] limma_3.42.2 fansi_0.5.0 pillar_1.6.2
## [145] fastmap_1.1.0 httr_1.4.2 survival_3.2-12
## [148] glue_1.4.2 zip_2.2.0 png_0.1-7
## [151] bit_4.0.4 class_7.3-19 stringi_1.7.3
## [154] sass_0.4.0 blob_1.2.2 latticeExtra_0.6-29
## [157] memoise_2.0.0 e1071_1.7-8 ape_5.5